Dephasing in coherently-split quasicondensates 
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We numerically model the evolution of a pair of coherently split quasicondensates. A truly 
one-dimensional case is assumed, so that the loss of the (initially high) coherence between the 
two quasicondensates is due to dephasing only, but not due to the violation of integrability 
and subsequent thermalization (which are excluded from the present model). We confirm the 
subexponential time evolution of the coherence between two quasicondensates oc exp[— (t/to) 2 ^ 3 ], 
experimentally observed by S. Hofferberth et. al, Nature 449, 324 (2007). The characteristic 
time to is found to scale as the square of the ratio of the linear density of a quasicondensate to 
its temperature, and we analyze the full distribution function of the interference contrast and the 
decay of the phase correlation. 
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I. INTRODUCTION 

Dephasing and dccohcrcncc arc phenomena at the 
heart of many-body physics which are deeply related 
to such fundamental problems as the crossover between 
quantum and classical dynamics of complex systems 
reversibility of physical processes iH or quantum infor- 
mation storage and processing |3j. To better under- 
stand dephasing and decoherence phenomena, we need 
systems, which are, on one hand, simple and theoretically 
tractable, but, on the other hand, available experimen- 
tally. In particular, ultracold, weakly-interacting bosonic 
atoms offer such an opportunity. In the present paper, 
we investigate numerically the time-dependent dephasing 
of ultracold atomic systems. 

A system of identical bosons confined to one- 
dimensional (ID) geometry is experimentally realizable 
with ultracold atoms trapped on atom chips Q or in an 
array of tight waveguides formed by a 2D optical lattice 
Q. The conditions of one-dimensionality are smallness 
of the temperature with respect to the energy quantum 
?iuj t of the radial (harmonic) Hamiltonian, 



(1) 



and smallness of the product of the 3D s-wave atomic 
scattering length a s and the mean linear density n^, 



ni D a s < 1. 



(2) 



A ID system of identical bosons interacting via contact 
pairwise potential is describable in the second quantiza- 
tion representation by the Hamiltonian 



n= dz 



f d ff + 9 ^mA, (3) 

2m oz oz 2 / 



where m is the mass of the boson, i[> = i>(z,t) is the 
bosonic annihilation field, and giD is the effective cou- 
pling strength of the ID contact interaction (in what fol- 
lows we assume repulsive interaction, giD > 0). The 



system described by this Hamiltonian supplemented by 
periodic boundary conditions is known to be fully inte- 
grable with the ground state properties and excitation 
spectrum in the thermodynamic limit given by the well- 
known Lieb-Liniger model @. 

Recalling that the coupling strength satisfies #id = 
2hoj I a B as long as a s -C l r = \/%/Jjrwx^) Q, we see that 
Eq. ([2]) requires smallness of the mean interaction energy 
per particle with respect to fkj T . Both Eqs. (fTJ) and @ 
mean small population of radially excited modes, either 
by temperature or interaction effects, respectively. 

In fact, the radial degrees of freedom are always excited 
virtually, with the excitation amplitude ~ niY>a B . This 
leads to the emergence of higher-orders in ip'ip (cubic 
etc.) terms in Hamiltonian, in addition to what is given 
by Eq. ©. These terms, corresponding to many-particle 
(three-particle etc.) effective elastic collisions, violate the 
integrability and lead to thermalization on the time scale, 
at longest, ~ l/[uj r (ni^ ) a 2 /l I ) 2 ] Q. Virtual radial mode 
excitations have been studied even earlier, in the context 
of soliton decay @ or quasi-lD (macroscopic) flow of a 
degenerate bosonic gas through a waveguide [l(|. How- 
ever, here we neglect this effect in order to study purely 
integrable dynamics. 

In what follows, wc consider uniform (nio = const 
for all z) and weakly interacting {mgm/h 2 <C rim) 
systems, tiid = (ip 1 !j(z,t)'ipj(z,t)) being the mean lin- 
ear density of bosons. For temperatures below T qc ~ 
(giD^nf^/m) 1 / 2 /&B [H|, a weakly-interacting system 
of bosons is in the quasicondensate state [12} , which 
means that the operator ip may be replaced by a clas- 
sical complex-valued field with a phase fluctuating along 
z (density fluctuations in the practically important long- 
wavelength range are suppressed via interactions) . In this 
regime, not only the phase coherence is mantained, but 
also the density-density correlation function at zero dis- 
tance approaches 1 (instead of 2, the value characteristic 
for a non-degenerate bosonic gas) [ll|. The stationary 
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two-point correlation function of a quasicondcnsate at a 
finite temperature is [l3[ 

^}(z,t)^(z',t))=n 1D cxp(-\z-z'\/X T ), (4) 

where the thermal phase-correlation length is 

A T = 2h 2 n 1D /(mk B T). (5) 

Quantum noise is dominant on length scales shorter than 
)] 1 /2 « At [3. Finite-temperature 
fluctuations in degenerate bosonic gases in highly 
anisotropic traps has been extensively studied both the- 
oretically [HI EBl and experimentally [13, E3 • 

One of the fundamental questions related to the Lieb- 
Liniger model is: how fast will two mutually decoupled 
quasicondcnsates with equal mean linear densities Uyd 
decohere, if their initial fluctuations are highly correlated 
(ipl(z, 0)ip2(z, 0)) « niD? Each of the quasicondensates 
has (upon tracing out the variables of another one) a 
finite temperature T. The measure of coherence will be 
the coherence factor 

y(t) = (Re{exp[i(<p 1 (z,t)-<p 2 (z,t)-6 ov (t))} }), (6) 

where (pj(z,t) are the local phase operators for the 
condensate labeled j. By 6 ov (t) we denote the con- 
stant overall phase of the system ( the phase associ- 
ated with the Goldstonc mode) defined by 9 ov (t) = 

2u r a s f* dt' Jq dz [|-0 2 (-Z, i')| 4 — l^i [z, t')\ 4 ] which was 
described by Lewenstein and You [ljf (with L the to- 
tal length of the condensate). From now on, we replace 
operators ipi, ip2 by complex random functions ipi, ip2, 
whose fluctuations account for thermal noise. 

We use "operational" definition © of the coherence 
factor on the two reasons: (1) density fluctuations are 
suppressed by interatomic repulsion and therefore Eq. 
(|6]) is practically equivalent to the more traditionbal def- 
inition *(t) = n^{4>2(z,t)4!i(z,t)e~ ieovW } (accounting 
for the correction to the overall phase); (2) in a time- 
of-flight interference experiment, the relative phase is di- 
rectly measured, while the density fluctuations are much 
harder to detect. We retain the symbol of the real part in 
Eq. ([6]) , because (i) the imaginary part of the expression 
in curly brackets becomes exactly zero only after averag- 
ing over an infinite ensemble of realization and is small 
but finite for real experimental data; (ii) we want to be 
consistent with the notation of Rcfs. (20. |2lj. 

The phase and density fluctuations in a quasiconden- 
sate are calculated from the harmonic approximation to 
the exact Hamiltonian In the harmonic approxi- 

mation, fluctuations are Gaussian with zero mean, and, 
hence, Eq. © reduces to 

*(t) = exp{~{[<p 1 {z,t)-(p 2 {z,t)-6 ov {t)} 2 )\. (7) 

In experiment [22} . a subexponential decay of coher- 
ence 

*(f) « cx V [-(t/t a ) a ] (8) 



has been detected, with the numerical value of a statis- 
tically consistent with the hypothesis 

a = 2/3. (9) 

Initially there was a single quasicondensate of 87 Rb atoms 
at the temperature T- ln between 82 nK and 175 nK. Then 
it was coherently split into two quasicondensates with 
the density n^o each (the values of nro were in the range 
from 20 fim^ 1 to 52 ^inr 1 ). The splitting was made 
as adiabatic as possible, so that the initial fluctuation 
in both quasicondensates just after their full separation 
were almost identical. However, the coherence factor 
then decayed rapidly as the hold time grew, according to 
Eqs. ([HIO, with to ~ 0.01 s, which is an order of magni- 
tude shorter than the expected thcrmalization scale for a 
non-degenerate system of that density Q. This obtained 
to was in a fair agreement with the theory developed by 
Burkov, Lukin, and Demler [2(J, which will be briefly de- 
scribed later, under the assumption that the temperature 
of the two quasicondensates after the splitting was ap- 
proximately equal to Tj n . However, the range of parame- 
ter variations (the radial trapping frequency was chosen 
to be either 2tt x 3.3 kHz or 2tt x 4.0 kHz) in Rcf. [H| 
is too narrow to reliably determine the dependence of to 
on n-iD, T, and cj r . 

II. THEORETICAL APPROACHES 

The two existing theoretical descriptions [2(| [2l[ of de- 
phasing in ID quasicondensates share the common basic 
model but differ in technical tools to solve it. In both 
cases two quasicondensates, describable by Eq. ©, arc 
separated by a potential barrier wide and high enough to 
make their tunnel coupling negligible and, hence, their 
time evolution after splitting fully distinct. These quasi- 
condensates are assumed to have the same mean atomic 
density and placed in two waveguides with the same ra- 
dial trapping frequency. The temperature is low enough 
{k^T < fi = 2?iuj T niT)a s ) to consider only phononic part 
of the elementary excitation spectrum. Although static 
Eqs. (|U [SJ hold also for k^T > /i, the dynamic theo- 
ries of Refs. [10, [2l[ rely on the phononic type of the 
elementary excitations spectrum. In the case k^T > fx 
thermally populated particle- like modes contribute to the 
system dynamics, and we expect therefore a deviation of 
the low of coherence factor decay from the theoretical 
predictions [13, [H|. 

The quantum noise is fully ignored, and excitations are 
represented by small-amplitude classical waves. Initially 
the fluctuations in the both quasicondensates are almost 
perfectly correlated (initial small interwell fluctuations 
serve as a seed noise, whose detailed properties are not 
"remembered" by the system in the long-time asymptotic 
regime and thus do not affect the final result). 

The theory by Burkov, Lukin, and Demler poj was 
based on derivation of a Langcvin-type equation with 
the random source term correlation and damping-term 
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properties described by a certain kernel obtained by 
re-summation of diverging diagrams describing the ex- 
change of quasiparticles between symmetric ($ + = {ipi + 
■02)/a/2) and antisymmetric ($4. = ($1 — V-'2)/v / 2) modes. 
Finally, Burkov, Lukin, and Demler obtained Eqs. ([8j [9]) 



with t = t BLD , 



BLD 



2MTrh)j,lC/{k B T) 2 , 



(10) 



where /C = ir^/hnm/(2u! r a s ) is the Luttingcr-liquid pa- 
rameter of the system. Initially, the symmetric mode is 
mostly populated, its initial temperature being T + | t= o ~ 
2T, and only very small number of excitations are present 
in the antisymmetric mode, because of slightly nonadia- 
batic splitting. At large t, the temperatures of two modes 
equalize, T + w T_ ps T. Thus, in a strict sense, Ref. [2(| 
reproduces Eqs. (J8j [9j in the long-time asymptotic limit 
only, although in experiment this subexponential coher- 
ence decay law was observed even for t < to- 

The Burkov, Lukin, and Demler theory has, although, 
a weak point: it predicts overdamping of modes with en- 
ergies larger than /i 2 /(k B TJC). The damping rate of such 
modes is of the order of or larger than their frequency. 
Probably, this result is associated with a possible tech- 
nical overestimation of the re-summed divergent series. 
It also may stem from the assumption of the purely lin- 
ear dispersion law for elementary excitations [20| that 
provides a large phase space available for products of a 
splitting of one phonon in the + mode to two phonons 
in the — mode. However, higher order corrections to the 
linear phononic dispersion law make this process energet- 
ically forbidden in ID via a small but unavoidable energy 
mismatch. Neglecting the latter fact may also results in 
obtaining too fast dynamics of the system. This moti- 
vated us to reconsider the problem and to put forward 
an alternative explanation of the subexponential dephas- 
ing. 

We refer a reader to our previous paper 



21| for the 



details of the calculations, which arc briefly summarized 
below. We considered motion of pairs of compact wave 
packets of phonons (with the localization size of the order 
of the carrier wavelength 2ir/k) in two random ID me- 
dia with relative fluctuations of local parameters (density 
and flow velocity) and follow their mutual dephasing, as- 
cribing the obtained dephasing rate to the elementary 
mode with the momentum hk. The source of mutual de- 
phasing is the dependence of the phonon frequency lu^ on 
the local density Sn and flow velocity SV fluctuations: 



dc 



dn 



ID 



5n + 6V \k\ 



(11) 



where c = y/ fi/m cx v/niD is the speed of sound. Then 
we assume that the contribution into the right-hand-side 
of Eq. (|lip comes only from fluctuations with the wave- 
length longer than 2ir/k (the influence of short-range fluc- 
tuations is averaged out). In other words, a propagating 
wave packet "sees" the fluctuations only on the length 
scales longer than its carrier wavelength. 



We calculate the statistical properties (the local cor- 
relation function at two different instants of time) of the 
differences of the fluctuations related to the 1st and 2nd 
quasicondensates, i.e., of 5n\ — Sri2 and 5V\ — 8V2 in the 
limit of asymptotically long time, when the symmetric 
and antisymmetric modes mostly equilibrate. The esti- 
mation for the dephasing rate for the two wave pack- 
ets with the momentum hk propagating in two parallel 
quasicondensates gives 



c\k\ 



dk' k B T 



8 J_ k 2ir /xniD 



1/2 



k B T 



fcl 3 / 2 . 



(12) 



Here we estimate <; ~ ^/5/ (8~7r) ~ 0.446. This estima- 
tion stems from our sharp-cutoff assumption. A different 
model, using some smooth function to eliminate the influ- 
ence of modes with momenta 3> hk, would give another 
value for However, later we shall see that <r ~ 0.446 is 
quite a reasonable value. 

Note that Eq. (Ti"2"|) does not predict overdamping of the 
modes with the energies close to k B T < \x: their damping 
rate is less than their frequency by a factor ~ /C. 

The shorter the wavelength of a mode, the faster this 
mode equilibrates. The integrated coherence factor is 
then 



ty(t) = exp <^ - 



mk B T 
2irh 2 niT> 



dk k~ 2 [1 - exp(-r fc t)] 



(13) 



Since Tk oc \k\ ' we obtain by integration ^(t) 
exp[-(t/*o) 2/3 ] but with to = t^ s [HI, 



3^2 



h 6 n 



ID 



m{k B Tf 



k mX 2 -. 
4~h~' 



(14) 



where k w 2.85 if we take q = 0.446. The fact that t^ s 
does not depend on c and, hence, on the atomic inter- 
action strength in the limit K, 1, seems to be deeply 
related to the independence of the thermal correlation 
length At [Eq. (0] on c in the static regime. 



III. NUMERICAL MODEL 

We directly simulated the time evolution of two 87 Rb 
quasicondensates by solving (by a split step spectral 
method [23j ) two Gross-Pitaevskii equations with ini- 
tial conditions chosen randomly corresponding to Bose- 
Einstein statistics of classical (thermal) excitations in 
phase and density waves. The local phase <f>j and density 
Uj, j = 1,2, values in the 1st and 2nd condensates are 
expressed through the respective values for the symmet- 
ric and antisymmetric modes (the local velocity is related 
to the phase as Vj = (K/m)d<j)j/dz). At t = they are 



01,2(2,0) 



<j>+{z,Q)±<j>-(z,0) 
T2 ' 



4 



/ ~\ Sn + (z, 0) ± 8n^(z, 0) . . 
n 1)2 (z,0) = n m + +v ^ (15) 

where the fluctuations for the + and - modes will be 
specified below. 

A. Initial conditions: The splitting process 

Exact values of these fluctuations depend on the details 
of the splitting process and its non-adiabaticity (for ex- 
ample, for the regime of the linear decrease of the tunnel 
coupling between two quasicondensates the population of 
phonons in the - mode is expressible via Besscl functions 
|24|). However, without exact knowledge of the splitting 
process, we can only give an estimation of the initial fluc- 
tuations. In the present paper, we approximate the initial 
populations of the symmetric and antisymmetric elemen- 
tary excitation modes by thermal distributions with the 
temperatures T + and T_ = t]T+, r\ -C 1, respectively. In 
the course of subsequent evolution, the temperatures of 
the + and - modes equalize and approach [(1 + r/)/2]T + . 

The particular choice of this parameter (r\ = 0.1 in 
our simulation) and the stability of our results against 
its variations are discussed in Sec. IV. B. 

The initial fluctuations appearing in Eq. (|15[) are given 

by 



' k 

0±(z,O) = -±—'£*B± S m(kz + £),(16) 

where Sk = \k\/ \J k 2 + Amfi/h 2 . Here the sum is taken 
over k = 2irMk/l, where Mk is a non-zero integer 
number running from — M max to M max . The phases 
arc uniformly distributed between and 2ir and 
is a positive number whose square has the ex- 
ponential probability distribution dP(\B^\ 2 )d\B^\ 2 = 

<|B±| 3 )-iexp(-|fl±| 2 /<|B±| 2 )). 

In the present paper, we completely neglect the quan- 
tum noise (zero-point oscillations of quantized local den- 
sity and phase), since taking it into account severely 
limits the time scale of reliable numerical integration 
of equations of motion of the system in the truncated 
Wigner approximation [25j . Since we take into account 
only classical (thermal) excitations, we have: 

(\b£\ 2 )= i kBT± (17) 

Here is the difference between our approach and that of 
Bistritzer and Altman [26| , who simulated the dephasing 
of two quasicondensates with only quantum (T_ = 0) 
fluctuations in the antisymmetric mode at t = 0. Us- 
ing pseudo-random numbers and £j^ 2 ' uniformly dis- 
tributed between and 1, we obtain the amplitude and 



phase of the excitation in a given mode: 

Bt = yJ(\Bt\ 2 )\l*tU <*=2< a . (18) 

B. Numerical method 

For obtaining reliable random numbers for the choice of 
initial states, we employ a Pseudo- Random Number gen- 
erator of Wichmann and Hill [27| , implemented here as 
a 4-fold combined multiplicative congruential generator. 
We use periodic boundary conditions with a simulation 
interval length which is of the order of the longitudinal 
extension of the condensate. 

Then the time-dependent Gross-Pitaevskii equation is 
solved by a 4th-order Fourier Split-Step method, the cor- 
responding time-dependent coherence factor is calculated 
and then averaged over a large enough sample of ran- 
dom realizations of this simulation. The resolution of the 
Fourier approximation is chosen such that the maximum 
energy of the quasiparticle (corresponding to the shortest 
resolvable wavelength) exceeds k-gT, thus preventing the 
overestimation of the intcrwcll coherence. 

We obtained reliable (stable and convergent) results 
in the parameter range spanning over half an order of 
magnitude in both the temperature T (from 30 to 100 
nK) and the mean densities Uid (from 30 to 100 /im -1 ) 
with and additional constraint (to be discussed below): 
the reliable results were obtained for the density-to- 
temperature ratio less than or of order of 1 nK -1 /im _I , 
further increase of this ratio resulted in the apprecia- 
ble dependence of the results on the grid size, which we 
could not eliminate within the reasonable range of the 
space and time steps. We also varied the effective ID 
coupling strength by increasing the radial trapping fre- 
quency froni27r x 3 kHz to 2tt x 9 kHz. Typical results 
are shown in Fig. [T] 



IV. NUMERICAL RESULTS 

A. Full distribution function of interference 

In Fig. [3J we illustrate the numerically obtained tem- 
poral evolution of the two-condensate joint Full Distri- 
bution Function (FDF), in the spirit of Ref. [Hj]. On 
polar plots shown in Fig. [2l the polar radius of a point 
gives the contrast C of the interference fringes integrated 
over the sampling length £ sa m, and the polar angle is the 
phase O of such an integrated interference pattern. The 
FDF (in arbitrary units) is shown as a false color density 
plot. 

Locally, the relative phase randomizes rapidly, we can 
see that from Fig. [2] (a). However, since the correla- 
tion length of the phase in each quasi-condensatc is At, 
the contrast does not decrease significantly as long as 
isam -Vr. As we increase the sampling length, the 
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FIG. 1: (Color online). Typical distribution of the relative 
phase 6(z,t) = <f>i(z, t) — <j>2(z,t) along the quasicondensate 
axis at several hold times between t=0 and t = 0.25 s. nm — 
50 Atm" 1 , T+ = 50 nK. 

local phases become more and more averaged out, and 
for the parameters of Fig. [2] (b) , the contrast C begins 
to decrease first, and then the distribution of 6 starts to 
spread. Remarkably, even at times as long as 0.25 s there 
are still many realizations yielding C close to 1. 

Note, that in our analysis we assumed equal mean 
atom numbers in the two quasicondensate. The splitting 
process, that always happens in finite time, causes fluc- 
tuations of the relative number difference between two 
wells and, hence, provides an additional mechanism for 
the global phase diffusion [2^, [3(j. Experimentally ob- 
served [22j high initial phase coherence signifies uncer- 
tainty of the interwell atom-number difference (otherwise 
the overall phase would be completely random). How- 
ever, the phase diffusion affects only the global phase 
and not the coordinate-dependent noise and correlation 
properties. The global phase can be eliminated during 
the elaboration of experimental data and is therefore not 
a major hurdle to experimental studies of dephasing. 



B. Evaluation of the coherence factor 
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FIG. 2: (Color online). Time evolution of the Full Distribu- 
tion Function shown as a false-color density plot (dark-blue: 
zero density; red: maximum density) derived from 3000 nu- 
merical runs. Hid = 50 fj,m~ , T+ = 50 nK. Increasing time 
from top to bottom, times between s and 0.25 s. Sampling 
length: L sam = 6 /im for left column, I/ sam = 50 fim for right 
column. Units on the axes are dimensionless. 



Evaluating our simulated data, we obtained a subex- 
ponential decay of the coherence factor consistent with 
Eqs. (JSJ in a wide range of parameters, see Fig. [3] 

Before proceeding further, we discuss our choice of the 
parameter 77, the initial ratio of the temperatures in the 
+ and — modes, and the sensitivity of our results to 
changes of this parameter. 

For all the simulations, presented beginning from Fig. 
21 we used r\ = 0.1. We choose this value because it 
gives the initial contrast that agrees with the experimen- 
tal data [iH [^(0) ~ 0.9 for a condensate of the length 



~ 100 /im]. Anyway, a description of the initial fluctu- 
ations beyond our model of thermal fluctuations at the 
temperature T_ = r]T + would require a detailed knowl- 
edge of the splitting process. 

We checked the sensitivity of our simulations to the 
choice of the parameter 77. A plot of In | In <3> | against \nt 
is shown in Fig. [3] We can observe that the slope of 
the sub-exponential decay is varying only slightly. This 
means that the subexponential decay with a — 2/3 is 
quite insensitive to initial conditions within the given 
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FIG. 3: (Color online). In | In 'S' | as a function of time for 
different values of 77 = ^F: 77 = 0.14 (red solid line), 0.1 (blue 
dashed line), 0.08 (cyan dash-dotted line), and 0.06 (green 
dotted line). Logarithmic scale for time. Straight solid line 
with crosses at the ends shows slope of 2/3. tiid = 50 /im" 1 
and T+ = 70 nK for all curves. 



range of 77. However, it remains unclear why this regime, 
which appears as asymptotic in both existing theories 
[20I l2lj , sets on after a very short time (few ms) . 

We now turn to the question of the decohcrcnce time, 
and how it scales with the different parameters of the Id 
system. In Fig. Ul we show numerically obtained values 
of to as a function of riirj. In Fig. [5j we show a plot 
of to against the ratio tiid/T + . We discern two different 
ranges of parameters. If niu/T + < 1 nK - /im -1 , we 
observed a sub-exponential decay of <F with the value of 
a from the interval between 0.66 and 0.69. The corre- 
sponding decay time is fitted by the formula 



6.4 



ft 3 "ip 
m{k^,T+)' 2 



1.6^ 



(19) 



T=T + 



Since in our calculations we took 77 = 0.1, and, hence, 
the temperature of both (+ and -) modes close to their 
equilibration is T sa 0.55 T+, Eq. (fTi?|) corresponds to 
Eq. (|T4"|) with k w 1.9, which agrees by the order of 
magnitude with k ~ 3 approximately evaluated in our 
theoretical model. Fig. 2] shows the fitting of to by Eq. 

For tiid/7+ <^ 1 nK /*m , our numerical method 
starts to fail. We are still able to fit the contrast decay 
by the formula (|5J) with a > 0.7, but the corresponding 
values of to do not obey Eq. (|19p anymore, and we have 
instead to ~ (m/h)X^ 67 b - 33 , where the length b depends 
on the grid size, thus indicating a numerical artefact. A 
possible explanation of such a behavior is that for small 
temperatures and high densities the interwell coherence 
persists for a long time, and the numerical error in sim- 
ulations accumulates faster than the actual "physical" 
dephasing occurs, thus yielding the dephasing time de- 
pendence on the grid size. This explanation, yet not 
fully corroborated, is at least in accordance with the ac- 
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FIG. 4: (Color online). Dephasing time to as a function 
of mean density niD for different values of temperature T+. 
Symbols: simulation results, thin line: Eq. (|19|) . T + — 70 nK 
(blue asterisks and cross), and 90 nK (red squares). 



celerated (a > 0.7) decay of simulated coherence in the 
problematic parameter range. The blue cross in Fig. @] 
corresponds to a parameter set from this regime. 

In the case where /sbT+ > to starts to deviate from 
Eq. (Unj). The reason is that the theory resulting in 
Eq. (fl"4"|) or (p~9|) is based on the assumption that only 
the phononic part of the Bogoliubov spectra is occupied. 
This case corresponds to the first two leftmost values of 
Fig. El 

In the range of applicability of Eq. (fT9|) , we checked the 
independence of our result on the interaction strength, 
see Fig. [6] for a plot of 1 F for different values of radial 
trapping. The only differences visible in this comparison 
result from statistical fluctuations. This independence of 
the time evolution of the coherence factor on the inter- 
action strength confirmes the validity of the theoretical 
model proposed in Ref. 21 [ . 



C. Correlation function 

Additionally, we investigated correlation properties of 
the interwell coherence, which are an important tool of 
understanding the quasicondensate properties 0, HH, 
|32| |. Knowing the interwell coherence autocorrelation 
function one can, under Gaussian fluctuation assump- 
tion, numerically construct distribution of the coherence 
factor (and the associated phase of the integrated inter- 
ference pattern) for any sampling length. In contrast to 
these previous works, we deal now with non-stationary 
correlation properties. We express the autocorrelation 
function of the interwell coherence via a new function 
<f>(z - z',t) as 



(ri(z)Mzm(z')Mz')) = n? D exp(-4>), 



(20) 



<f> being determined mainly by phase fluctuations (den- 
sity fluctuations are suppressed by inter-atomic repulsion 
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FIG. 5: (Color online). Dephasing time to as a function of 
the ratio 71id/7+. Solid blue line: Equation (I19|) . Symbols: 
simulation results, up triangles: T+ = 50 nK, square: T+ = 60 
nK, open circle: T+ = 65 nK, down triangles: T+ = 70 nK, 
cross: T+ = 75 nK, filled circle: T+ = 80 nK, diamonds: 
T+ = 90 nK. Inset: the ratio of numerically obtained values 
of to to the values of to given by the fitting Eq. ()19|) . The 
deviation of this ratio by 15% from unity (horizontal dashed 
line) for the leftmost point is due to temperature high enough 
to significantly populate particle-like states and thus violate 
the assumption of the phononic dispersion law underlying Eq. 
(|19|) . The next point shows the same tendency, but to less 
extent. 



in the phononic regime) . Analogously to Eq. ([7]) , we may 
write 

$ = -^ 1 (z,t)-^ 2 (z,t)-v 1 (z',t) + V2(z',t)} 2 )- (21) 

We calculate the right-hand-side of Eq. ([21]) by gener- 
alizing the formula obtained by Whitlock and Bouchoule 
[3l| . where we set the interwell tunnel coupling to zero 
and assume that, in the transient regime, mean popu- 
lation of each antisymmetric mode is given by its own 
temperature Tk(t), which depends on the mode momen- 
tum fik and evolves in time approximately as [2lj 



(22) 



[cf. Eq. (fl~3]) . where 77 = <gC 1 is neglected]. Thus we 
write 



' dk mk * Tk{t \ l~cosk\z- z'\). (23) 

h 71 ID 



We expand cos k\z—z'\ in series in powers of its argument, 
perform integration over k for each term expressing the 
integrals via the gamma-function T(s) = J Q dx e~ x x s ~ 1 
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0.2 




°* no 



0.1 0.2 0.3 0.4 

t(s) 



0.5 



FIG. 6: (Color online), ^(t) for several values of radial 
trapping frequency ui r , nm = 50 (im _1 and T+ = 70 nK. 
oj r = 2ir x 3 kHz (blue diamonds), 2ir x 6 kHz (green aster- 
isks), and 2tt x 9 kHz (red squares). 



and re-assemble the obtained terms into hypregeometric 
series 



iFrn ipl 1 ■ ■ • i ; Ci , . . . , C m , x) 1 ~h 



bib 2 



bn X 



C1C2 ...Cm 1! 

&l(frl + l)6 2 (fr 2 + l) ■■■ 6n(bn + l) X 2 

ci(ci + l)c 2 (c 2 + l) ...c m (c m + l) 2! 
Finally, we arrive at the following expression: 



2\z-z'\ 
At 



4(1-7?) ■ 

1 -3^r+^ £a(l *-* I/a) 



(24) 



where At is to be evaluated at the equilibrated tempera- 
ture T = Too = {T+ +T_)/2. The new auxiliary function 
S depends on z — \z — z'\/a only, where 



a 



2+2 \ 1/3 



(25) 



mniD 

S can be expressed via hypergeometrc series as follows: 

1 5 11 1 5 7 4 



Z 
-A 



=m - 1 f^lK pf 1 5 n 1 5 7 4 . 4 

2 V 6 12 12 2 6 6 3 729 



~z e ) - 



Z / 1 

24 4F H2' 

(-;) 



1 3 



5 5 7 4 3 5 



72!) 



S 6 ) + 



r ( — i x 

3/ 

1 J? ( 1 1 7 • 1 1 2 1 

V 6 12 12 6 2 3 6 



(26) 



729 / J J 



$ changes significantly on a time scale to. At i = 
we have $ = [4r?/(l + rfflX^lz - z'| < \z - z'\/X T . In 
the opposite limit t 7j we have a — > oo, 5 — > for 
any finite coordinate difference, S — > 0, and therefore 
$ = 2|z- z'|/A T . 

An important advantage of the autocorrelation func- 
tion (|2"Tj)) is its independence of the global phase diffusion 
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FIG. 7: (Color online). Function $ (dimensionless) that de- 
fines the autocorrelation function of the interwell coherence 
Eq. Q20p vs. the co-ordinate difference. Special symbols: nu- 
merical simulation results. Lines: fitting with ? = 1.0. The 
quasicondensate parameters are nm = 40 /Ltm , T+ = 30 
nK. The time elapsed after the end of coherent splitting pro- 
cess is t = 0.1 s (asterisks, dashed line), 0.2 s (diamonds, solid 
line), and 0.5 s (crosses, dot-dashed line). 



[19t l29t l3Cj ■ Therefore the theoretical predictions for the 
autocorrelation function (|20p admit direct comparison to 
experiment, unlike the contrast "J, which needs a cor- 
rection to the phase diffusion. If the latter effect can 
not be neglected, the direct comparison of experiment to 
theory in terms of the contrast "J is hindered, because un- 
ambiguous reconstruction of the global phase shift from 
the measurements in each experimental run is practically 
very difficult. 

In Fig. [7]we display the simulated values of $ and their 
fitting by Eqs. (|24l - [26)) at various times t. The only 
fitting parameter is <; whose value is found to be £ « 1.0, 
which is in a fair agreement with the value <r « 0.74 that 
corresponds to the numerical prcfactor in Eq. (|19|) and 
is not too far from the estimated value <; = 0.446. 



D. Comparison to experiment 



In the experiment [22( the temperature T n before split- 
ting was known, and the temperature Tf after splitting 
and equilibration between the + and - modes was es- 
timated using the theory of Ref. [2(|. Now we repeat 
the same procedure using our Eq. (|T9|) . Recalling that 
i| « 1, we estimate Tf rj 0.5 T + , where T + is obtained 
from Eq. (|19p where to is given by the experiment for 
certain T- m and nm and present the results in Fig. [3] 
Surprisingly, the values of Tf estimated from two models 
(20T ] and [2l[ are quite close to each other and are simi- 
lar to T n . The radial frequency range employed in (22j 
(from 27T x 3.3 kHz to 27r x 4 kHz) is too narrow to allow 



4- 

H 



200 



150 



100 



50 



► 

> 



4. 



50 



100 



150 

(nK) 



200 



FIG. 8: (Color online). Final temperature after equilibrat- 
ing in experiment extra pola ted by Eq. (|19[) (open triangles) 
and by theory of Ref. [2Q| (filled triangles). At radial fre- 
quency 2-7T x 3.3 kHz (smaller symbols): background density 
?iid = 20 (left triangles), 34 (up triangles), and 52 fim -1 
(right triangles). At radial frequency 2n x 4.0 kHz (larger 
symbols): background density tiid = 22 (left triangles), 37 
(up triangles), and 51 ^raT 1 (right triangles). 



unambiguous determination of the dependence of to on 



Relation between T+ and Tin under adiabatic splitting 
conditions 



Obtaining the relation between T+ and T n is a subtle 
matter, not covered by the existing theories. Here we 
propose a way to estimate the ratio T + /T; n from sim- 
ple considerations. The adiabaticity of the splitting pro- 
cess means that the entropy (but not the energy) is con- 
served. In other words, populations of the momentum 
states should not change in the course of splitting. Just 
after the splitting most of the noise occurs in the sym- 
metric mode. However, the speed of sound in the + mode 
changes, as we show below. 

Indeed, the speed of sound c is to be calculated from 
the formula 



47r?i 2 a, 



dx I dy\ty±(x,y)\ 4 n to t, 



(27) 



where "^(cc,?/) is the wave function (normalized to 1) 
of the ground state of transverse trapping Hamiltonian. 
To lowest order, we neglect the influence of the atomic 
interactions on the transverse profile of the quasiconden- 
sate (ToT |. in agreement with our assumption of a truly 
ID regime. Eq. (|2T|) applies for both single- and double- 
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well potential shape, n to t = 2«id being the total linear 
density of atoms in the system. 

The fundamental frequency of the transverse trapping 
potential is u r before splitting. Therefore initially 

*|f (x, y) = (V^kr 1 exp[-(x 2 + y 2 )/{2l*)]. (28) 

The splitting is designed so that in the end the funda- 
mental frequency of the radial oscillations in each of the 
two waveguides is again oj t . Therefore, 

*i s (z, y) « [*1(x -w,y) + ¥l(x + w, y)]/V2 (29) 



is to be substituted into Eq. (|27|) to determine the speed 
of sound after splitting. The separation 2w between the 
two waveguides is so large that the overlap between the 
wave functions localized near x = ±ui, y = is negligi- 
ble, as it must be in the zero-tunneling case. From Eqs. 
(|2"T1 - |2"9"|) one can easily see that the speed of sound drops 
after splitting by v2- 

As we discussed before, adiabatic splitting concerves 
the mode populations. If the temperature is comparable 
to the chemical potential, phonon-like modes (with the 
energy linearly proportional to the speed of sound) dom- 
inate. For each momentum hk the ratio hc\k\/T should 
be the same before and after the splitting (we denote the 
respective values of the speed of sound by c ln and c as ). 
Therefore, we expect 



T in 



1 

71' 



(30) 



However, T + can be affected by non-adiabatic effects, 
which are far beyond our ID approach, such as heat- 
ing via creation of vortices in the course of splitting and 
their subsequent decay [33| or heating by technical noise. 



Therefore it is extremely difficult to reliably predict a 
priori the ratio T + /T- ln and, moreover, Tf/T; n . 



V. CONCLUSION 

To conclude, we simulated numerically the evolution of 
two coherently-split quasicondensates in a broad range 
of experimentally feasible densities, temperatures, and 
effective ID coupling strengths (radial trapping frequen- 
cies). We reproduced the subcxponential decay Eq. dHJ 
of the interwell coherence with a very close to observed 
[lU and theoretically predicted [13, EH value 2/3. Our 
characteristic dephasing time to varies quadratically with 
the ratio of the linear density to the temperature (as long 
as fcsT is smaller than the chemical potential and thus 
only the phonon part of the elementary excitation spec- 
trum is thermally excited) and does not depend on the 
ID coupling strength. The latter fact has its static coun- 
terpart - the independence of thermal phase-coherence 
leng th At of a quasicondensate on the speed of sound c 
[13j . see Eq. ([5]). In other words, the typical dephas- 
ing time scales as to ~ mX^/fr, in accordance with our 
model (2lj . This conclusion is corroborated by analysis 
of autocorrelation of the interwell coherence during the 
dephasing process. Although comparison to the exper- 
iment [22( does not show the clear preference of one of 
the two theories [I(|[2l| over another, our numerical sim- 
ulations in a broad range of parameters unambiguously 
support our theory [2lj . 
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